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ABSTRACT 


In  the  operation  of  CMGs  there  exists  a  concept  called  “back  drive,”  which  represents  a 
case  where  the  coupling  effects  of  the  angular  velocity  of  the  body  and  the  angular 
momentum  of  the  CMG  overwhelm  the  input  torque  and  result  in  a  lack  of  control.  This 
effect  is  known  but  not  well  documented  or  studied  in  the  literature. 

Starting  from  first  principles,  this  thesis  derives  the  full  nonlinear  dynamical 
equations  for  CMGs.  These  equations  contain  significantly  more  terms  than  are  found  in 
the  literature.  As  a  means  to  understand  the  implications  of  these  terms,  a  reduced  order 
model  is  derived.  The  full  and  reduced  models  are  then  validated  by  means  of  extensive 
simulations.  Finally,  experimental  verification  of  the  models  confirms  the  finding  that 
the  reduced  order  model  provides  a  reasonably  high  fidelity  for  dynamics. 
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I. 


INTRODUCTION 


In  the  study  of  Control  Moment  Gyroscopes  (CMGs)  there  exists  a  concept 
called  “back  drive,”  which  represents  a  case  where  the  coupling  effects  of  the  angular 
velocity  of  the  body  and  the  angular  momentum  of  the  CMG  overwhelm  the  input  torque 
and  result  in  a  lack  of  control.  This  effect  is  known  in  industry  but  not  well  documented 
or  studied  in  the  literature. 

The  standard  model  of  a  spacecraft  with  Control  Moment  Gyroscopes  used  by 
most  textbooks  and  papers  assumes  that  the  control  variable  of  the  CMG  is  the  CMG 
gimbal  rate.  In  this  model  the  complete  angular  momentum  (h)  of  a  single  CMG  is 

assumed  to  be  the  moment  of  inertia  of  the  rotor  multiplied  by  the  angular  velocity  of  the 
rotation  about  the  spin  axis,  shown  in  Equation  (1)  where  J r  is  the  rotor  inertia,  Q  is  the 

rotor  angular  velocity  and  g ,  is  the  spin  axis  of  the  CMG,  as  shown  in  Figure  1 . 

h  =  Jrpgt  (1) 


1 


The  angular  momentum  of  a  system  of  n-CMGs  in  the  body  frame,  which  will  be 
fully  derived  in  Section  IV  but  is  useful  now,  is  then  defined  as  Equation  (2)  where 
h0  =  J  O  ,  8  is  the  gimbal  angle,  f3  is  the  pyramid  skew  angle  and  <f>  is  the  rotation  angle 
due  to  the  placement  of  the  CMGs  as  shown  in  Figure  2. 


h*=2>0 


c8ic(j)i-%8ic[:h(l)j 

c8rS(/)i+s8icf:Jc(/)i 

sSfi/3 


(2) 


whereof  =  cos  (£),  sJssin^),  ctj)  =  cos(^)  ,s(j)  =  sin(^)  ,c/3  =  cos(/?) ,  and 
s/3  =  sin(/?) . 


A  typical  n=4  system  is  shown  in  Figure  2  as  an  example  of  the  geometry.  The 
angles  [3 ,  (f)  and  8  are  labeled  on  representative  CMGs  but  are  not  labeled  where  their 

inclusion  would  reduce  the  clarity  of  the  figure.  The  axes  labeled  bvb2,  and  lx  constitute 
the  satellite  reference  frame. 
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Figure  2.  A  typical  n=4  CMG  Pyramid.  Both  <f)  (offset)  and 

1 3  (skew)  angles  are  fixed,  determined  by  the  mounting  of  the  CMGs 
with  respect  to  the  spacecraft  reference  frame  ( f\ ,  b2 ,  b3 ). 


With  the  angular  momentum  defined  as  in  Equation  (2),  for  any  given  CMG  h0 , 


P ,  and  <j>  are  all  constant  in  time  and  S  is  the  only  time  varying  component.  Therefore 


the  derivative  can  be  expressed  as: 


n 

f 

C(fcr/i(-S(f  c/A^( 

\ 

1-  =  A.Z— 

°t Idd 

V 

sSfiP 

) 

(3) 


If  we  define: 


3 


Then: 


a,  =  h0 


es 


cdy^-sSyfhcfi 
cdls(f)i+sdic/:jc(f)i 
sSjS/3 


J  ) 


h  =  K  •••  a„; 


8 


8. 


=  A5 


(4) 


(5) 


At  a  system  level,  the  rotational  equation  of  motion  of  a  spacecraft  with  a  CMG 
system  is  given  as  Equation  (6)  where  hn,v  is  the  angular  momentum  of  the  entire 


spacecraft-CMG  system. 


=  h  +  toxh 

sys  sys 


(6) 


Because  angular  momentum  is  conserved  we  can  write  the  angular  momentum  of 
the  system  as  the  sum  of  the  angular  momentum  of  the  spacecraft  body  and  the  angular 
momentum  of  the  CMG  system  as  shown  in  Equation  (7). 

h  sys  =  h*  +  h CMG  =  J*01  ^CMG  O 


where  3b  is  the  inertia  matrix  of  the  spacecraft  and  to  is  the  angular  velocity  of  the 

spacecraft  expressed  in  the  body  frame  (b  in  Figure  2).  Combining  Equations  (6)  and  (7) 
results  in: 

text  =  J b^h'N  +  ^CMG  +  X  (' J60)/’/JV  +  hCMG  )  (8) 

By  defining  the  internal  control  torque  generated  by  the  CMG  system  as: 

U  =  -(hcMG  +  0)WArXhCMG)  (9) 

We  obtain: 

4-  |  __  t  'b/N  .  b/N  x  ,-b/N  f\\ 

text  +  u  =  J  b(j)  +(0  X  JA(0  (10) 

Combining  Equations  (9)  and  (5)  the  gimbal  rate  can  be  obtained  as 

8  =  A_1(-u-©6WxhCMG)  (11) 
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Equation  (1 1)  is  only  valid  for  n  =  3  ,  therefore  a  more  general  approach  uses  the 
pseudo-inverse  shown  in  Equation  (12)  to  obtain  the  commanded  gimbal  rate. 

5  =  (Ar(AAr)  1)(-u-to6/A,xhCMG) 


This  results  in  the  state  space  equations  in  Equation  (13). 

q=|(®*q) 

«>  =  J*1  ((‘ext  +  U  )  -  X  J )) 

^  =  (Ar  (AAr)  1)(-u-ojb/N  xhCMG) 


Where  to*  is  defined  as: 


to  = 


0 

-co0 

cox 

-CO 3 

0 

Q)x 

co2 

a>2 

-CQy 

0 

co2 

-cox 

-®2 

-CU, 

0 

(12) 


(13) 


(14) 


This  7+n  state  space  model  is  sufficient  for  many  applications,  but  using  this 
technique  serves  to  mask  the  inner  dynamics  of  the  CMG.  Figure  3  shows  many  of  the 
intricacies  of  the  inner  dynamics  of  a  CMG,  which  are  neglected  in  the  analysis  shown 
above. 


GM 
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'iV«:ofdx 

H 

h- 

7] 
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Figure  3.  CMG  Inner  Loop.  After  Leonard.1 
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In  the  physical  world,  the  control  variable  for  the  CMG  is  the  motor  current  of  the 
CMG  gimbal  motor.  Approaching  the  problem  from  this  perspective  increases  the 
complexity  significantly  and  the  relationship  between  motor  current  and  gimbal  torque 
can  be  considered  fairly  well  known.  Therefore  a  better  analysis  of  the  spacecraft/CMG 
system  can  be  conducted  using  the  gimbal  torque  as  the  control  variable.  To  accomplish 
this,  the  system  must  be  analyzed  using  the  detailed  dynamics  of  the  CMG-spacecraft 
system. 

This  thesis  includes  a  summary  of  past  and  current  literature  on  the  dynamics  of 
CMG  clusters  with  emphasis  on  derivations  of  the  equations  of  motion  and  the 
implications  of  the  gimbal  torque  equations,  followed  by  a  detailed  derivation  of  the 
dynamics.  The  full  equations  of  motion  are  modeled  and,  using  simulation  and 
experimental  data,  the  physical  implications  of  the  equations  are  determined  and 
discussed. 


6 


Section  Notes 


1  Barry  S.  Leonard,  Control  Moment  Gyro  Systems,  Vol.  AE3818  Lecture  Notes  (Monterey,  CA:  , 

2002). 
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II.  LITERATURE  REVIEW 


A.  MOTIVATION 

Much  research  has  been  conducted  on  Control  Moment  Gyroscope  controlled 
satellites  since  they  were  first  developed  in  the  1950s.1-13  However,  limited  published 
research  can  be  found  on  the  detailed  derivation  of  the  full  dynamics  of  CMGs  and  of  the 
derivations  that  exist,  none  could  be  found  that  discuss  the  physical  implications. 
Searches  of  several  databases,  including  the  American  Institute  of  Aeronautics  and 
Astronautics  (AIAA)  and  the  National  Aeronautics  and  Space  Administration  (NASA) 
document  server,  revealed  a  series  of  trends  in  the  publically  available  published 
research. 

Early  research  in  the  1960s  and  early  1970s  focused  on  the  basic  issues  of  using 
CMGs  in  space — many  specifically  focusing  on  SKYLAB.2’3’4’5  This  research 
occasionally  included  a  detailed  derivation  of  the  dynamics  of  a  CMG/spacecraft  system 
but  the  developed  equations  do  not  appear  to  drive  CMG  design.6  Once  the  research  on 
the  SKYLAB  system  was  completed,  there  appears  to  be  a  gap  in  publically  available 
research  on  the  dynamics  of  CMGs  throughout  the  remainder  of  the  1970s  and  1980s. 

There  are  a  scattering  of  papers  in  the  1970s  and  1980s,  primarily  focused  on 
using  CMGs  for  fine  pointing  of  space  telescopes.7  These  papers  are  largely  focused  on 
hardware  design  and  steering  law  refinements.8’9’10  The  dynamics  of  the  CMG  were  not 
a  primary  focus  of  this  research  and  the  published  hardware  design  research  does  not 
emphasize  gimbal  torque  requirements. 

The  1990s  and  early  2000s  seemed  to  see  a  resurgence  of  CMG  dynamics 
research,  much  of  it  focused  on  either  precision  attitude  control  of  spacecraft  or  in  the 
relatively  new  field  of  using  CMGs  in  robotics.11’12’13’14  A  great  deal  of  the  robotics 
research  focuses  on  scissored  pair  CMG  systems,  which  involve  some  unique  dynamics 
of  their  own.15 

An  overview  of  the  published  research  is  included. 
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B. 


EARLY  LITERATURE 


CMG  analysis  and  research  in  the  1960s  varied  widely.  However,  several  papers 
include  a  derivation  of  the  dynamics  of  a  CMG  system.  Havill  and  Rate liff(  1964) 16 
derive  the  equations  of  motion  of  an  ideal  single  degree  of  freedom  CMG  in  their  paper 
regarding  twin-gyro  attitude  control  systems  but  as  their  proposed  system  cancels  the 
motor  torque  by  design,  they  do  not  go  into  depth  on  the  topic.  The  equation  they 
develop  for  gimbal  torque  is  shown  in  Equation  (15). 

CmG3  =  (^g3  +  )[®3  +  _  G^) 

This  equation  neglects  terms  created  by  the  body  rate  coupling  but  includes  the  rest  of 
what  appear  to  be  the  commonly  accepted  terms. 

Liska  and  Jacot  (1966)17  derive  the  dynamic  equations  for  a  linearized  simple 
CMG  in  their  paper  on  the  use  of  CMGs  and  reaction  jets  for  manned  space  station 
attitude  control.  They  use  a  simplified  derivation  to  find  the  gimbal  torque  shown  in 
Equation  (16). 

Cmg,  ~ DS  +  J rQ<x>g  =  ( Jg}  +Jr,) $  OG 

This  equation  includes  a  drag  component  but  does  not  include  the  body  rate  cross  product 
tenn  and  the  input  from  the  spacecraft  acceleration.  Liska  and  Jacot  state  that  they  are 
neglecting  “as  insignificant  those  torques  caused  by  rotation  of  the  small  gimbal  mass  at 
vehicle  rates.” 

Schindelin  (1968) 18  did  not  explicitly  derive  an  equation  for  gimbal  torque  in  his 
paper  On  Space  Vehicle  Attitude  Stabilization  by  Passive  Control  Moment  Gyros  but  did 
derive  the  equations  of  motion  of  the  CMG  about  the  gimbal  axis  which,  due  to  the 
passive  nature  of  the  system,  was  set  to  zero.  Shown  in  Equation  (17)  the  equation  is 
effectively  Equation  (16)  without  specifically  including  the  drag  terms. 

0  =  C*7  +  J  \d -\- J  Clco  (17) 

V  ft  ri )  n  si  v  ' 
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These  papers  are  the  only  ones  of  the  era  that  could  be  found  that  included  any 
discussion  of  the  gimbal  torque  motor  and  all  derived  the  gimbal  torque  equation  as  a 
tangent  to  their  major  purpose.  None  were  found  that  specified  the  assumptions  made  to 
develop  their  versions  of  the  equations. 


C.  MODERN  LITERATURE 


Published  papers  in  the  1970s  and  1980s  focused  primarily  on  steering  and 
control  law  improvements.  A  number  of  NASA  papers  were  also  published  that 
contained  descriptions  of  the  hardware  for  specific  CMG  systems,  although  nothing  was 
found  that  detailed  decision  design  drivers.19’20 


During  the  1990s,  most  research  continued  to  focus  on  steering  law 
improvements,  with  an  emphasis  on  time  optimal  reorientations,  including  research  that 
proved  that  an  Eigenaxis  maneuver  was  not  time  optimal  in  all  cases.21  The  late  1990s 
and  early  2000s  also  resulted  in  a  number  of  papers  on  the  possibility  of  miniaturizing 
CMGs  for  use  in  small  spacecraft.22’23’24  Notably,  Busseuil,  et  al.  published  a  paper  in 
1998,  which  actually  discusses  gimbal  motor  sizing  for  miniaturized  CMGs  25  However, 
they  begin  with  the  assumptions  in  Equation  (18)  and  do  not  explicitly  address  the 
derivation  of  the  equations  of  motion  or  what  terms  they  are  neglecting. 


tCMG,  ~{jg,+  ) ^  +  Jri ^ ®g2 

C'VG,  =  ’A,  ^  ( 'A,  +  & ) 


(18) 


From  Equation  (18)  they  derive  the  required  maximum  torque  required  from  the  gimbal 
motor  as  shown  in  Equation  (19). 


. MAX 
lCMG2 


+J, 


E 


MAX 


+  J  „  fi  CO 


MAX 

~&2 


(19) 


This  is  the  only  paper  that  was  found  that  explicitly  states  an  equation  for  sizing  a  CMG 
gimbal  motor. 

Many  papers  of  this  era  use  Equation  (18)  as  a  baseline  assumption,  without 
addressing  the  simplifications  made.  In  this  thesis  these  equations  will  be  known  as  the 
Reduced  Model. 
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D.  DYNAMICS  TEXTBOOKS 


The  majority  of  current  dynamics  textbooks  include  a  fairly  straightforward 
derivation  and  discussion  of  the  use  of  CMGs  as  a  momentum  management  system  but 
generally  fail  to  delve  deeply  into  the  internal  dynamics  of  the  CMG  system,  focusing, 
reasonably,  on  the  issues  of  singularity  avoidance.  In  Analytical  Mechanics  of  Aerospace 
Systems,  Schaub  and  Junkins,  26  include  a  derivation  of  the  rotational  equations  of  motion 
for  a  system  with  a  single  Variable  Speed  CMG  (VSCMG),  then  expand  it  to  an  N- 
VSCMGs  system.  This  derivation  is  in  support  of  an  argument  for  using  Variable  Speed 
Control  Moment  Gyroscopes  for  singularity  avoidance.  They  derive  the  resultant  torque 
of  movement  of  the  CMG  rotor  as  Equation  (20). 

>V  =S2Jri  (G  +  g|7G  +  (iVw) 

('ffco,  +  J r g i  ho  +  oyo)^  (y-f2  —  J )  +  J rCl{d  +  ) —  2JrSco2  j  (20) 

+s\jri  (k[oo  +  d)  +  {(!Yo1{  ji.  -Jri 

As  the  CMG  gimbal  motor  operates  only  in  the  gimbal  axis,  the  motor  torque  is 
shown  in  Equation  (21). 

tg  ={Jrt  [gld)  +  S)  +  (a) \co2[jri  -  J,.  ))-/,.  Q^)  (21) 

The  derivation  included  in  this  document  confirms  the  equations  that  Schaub  and 
Junkins  developed.  The  approaches  used  in  the  two  derivations  are  sufficiently  different 
that  a  high  degree  of  confidence  in  the  resultant  equations  is  justified. 
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III.  DEVELOPMENT  OF  HIGH-FIDELITY  MATHEMATICAL 
MODELS  FOR  A  SINGLE  CMG  SYSTEM 


A.  RIGID  BODY  DYNAMICS 


For  a  general  system,  the  rigid  body  state  equations  can  be  derived  as  shown  in 
Equation  (22). 27 


qx  =  —(^q4col  —  ( q30)2  —  q20)3 )  j 

q2  —  4®2  —  —  *?3®l)) 

=  '^'(^4(y3  ~  (#2®1  _  ^1^2)) 

04  =-^(?1«1+M+M) 


(22) 


ft), 

d>2 

OJ3 


•/  t  )  ^2^3 


(t2  +  (>^1  0/3  ) 

J  2 

—  —  (t3  +  (J2  —  Jx )  ru,®, 

»y  -> 


) 

) 

) 


This  system  assumes  that  the  applied  torque  ( 7) )  is  external  to  the  system.  In  the  case  of 
CMG  controlled  spacecraft,  this  assumption  is  not  valid  due  to  the  complex  coupling 
effects  of  the  system. 


B.  SYSTEM  DEFINITION 

A  CMG  controlled  spacecraft  system  has  four  frames  of  reference:  the  inertial 
frame,  the  body  frame,  the  gimbal  frame  and  the  rotor  frame.  Figure  4  shows  a  diagram 
of  a  CMG  with  the  gimbal  and  rotor  frames  labeled.  By  definition,  r\  and  are 
collinear. 
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Figure  4.  CMG  diagram  showing  the  gimbal  and  rotor  frames  of  reference. 
Figure  5  shows  a  system  diagram  with  all  the  frames  of  reference  labeled. 
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Figure  5. 


System  diagram  with  the  gimbal,  rotor,  body  and  inertial  frames  of  reference. 


For  the  following  derivation,  subscripts  identify  which  frame  and,  when 
appropriate,  which  axis  the  component  is  reference  to.  Superscripts  for  angular  velocity 
note  the  rotation  frames  of  reference,  for  instance  cbh'N  indicates  the  rotation  of  the  body 
with  respect  to  the  inertial  frame.  A  single  superscript  used  with  a  derivative  term 


indicates  the  frame  in  which  the  derivative  is  taken,  for  instance 


</(©WJV)V 


dt 


indicates 


the  derivative  of  the  rotation  of  the  body  with  respect  to  the  inertial  frame,  taken  in  the 
gimbal  frame.  Each  frame  is  rotating  with  respect  to  its  mounting  frame. 

Therefore,  the  angular  momentum  of  the  system  can  be  written  as  Equation  (23). 

K  =  K+hg+K  (23) 
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By  definition  any  change  in  the  angular  momentum  of  the  system  is  the  result  of 
external  torque  applied  to  the  system,  shown  in  Equation  (24). 


t 


ext 


dh,0, 

dt 


(24) 


If  the  external  torque  is  known  we  can  solve  for  d>  for  the  system,  which  can  then 
be  used  to  solve  the  system  dynamics. 

With  the  assumption  of  rigid  attachment,  we  also  assume  that  the  rotation  of  the 
CMG  gimbal  and  rotor  are  perfectly  aligned  and  can  be  defined  as: 

cf'8=Orl  cos,b=Sg 3  (25) 


The  inertia  dyads  for  each  component  of  the  system  are  defined  as  shown  in 
Equation  (26). 
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4. 

The  angular  momentum  of  each  component  is  therefore: 

hh=Jh-<»b,N  hg=Jg(os,N  hr=Jr-orW  (27) 


From  the  addition  theorem  of  angular  velocity,  we  can  show  that  the  angular 
velocity  of  the  gimbal  and  the  rotor  in  the  inertial  frame  are  equal  to  the  sum  of  their 
angular  velocities  in  each  applicable  frame: 


to 


g/N 


=  (O 


b/N 


+  (0 


gib. 


(O 


■IN 


=  (OhlN  +<o8,b  +(Or,s 


(28) 


For  the  purposes  of  the  initial  single  CMG  derivation,  we  assume  that  the  body 
frame  and  the  gimbal  frame  are  initially  aligned  and  that  the  only  variation  between  the 
two  frames  takes  place  about  the  gimbal  axis.  Therefore,  the  Direction  Cosine  Matrix 
between  the  gimbal  and  the  body  frame  can  be  defined  as: 
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(29) 


cos(c>)  sin  ( £)  0 

DCMg/A=Cg/6=  -sin(<S)  cos(S)  0 

0  0  1 

We  then  define  then  angular  velocity  of  the  body  frame  in  the  inertial  frame  as: 

cox 

a>2  (30) 

co^ 

to  can  then  be  transformed  to  the  gimbal  frame  as  shown  in  Equation  (31).  <obg/N  is 
defined  as  the  angular  velocity  of  the  spacecraft  in  the  body  frame  translated  into  the 
gimbal  frame.  This  step  is  necessary  for  future  calculations  in  the  gimbal  frame. 

cos(cS)  sin(cS)  0  cox  cox  cos(cS)  +  (y2sin(cS)  Q)g] 

0)hg,N  =  -sin(c7)  cos(J)  0  co2  =  -cox  sin(fi')  +  (u2  cos(J)  g=  cog  g  (31) 

0  0  I  ox  ox  (,)' 

With  these  definitions  in  place  we  can  separate  the  angular  momentum  of  each 
component  and  take  the  derivative  in  the  inertial  frame  in  order  to  determine  the  overall 
system  6)  the  inertial  frame. 


C.  EQUATIONS  OF  MOTION  FOR  THE  BODY 

As  defined  in  Equation  (27),  the  angular  momentum  of  the  body  is: 

h„=Jbo>h,N  (32) 

where  all  terms  are  defined  in  the  body  frame.  Using  the  transport  theorem,  the 
derivative  of  Hh  in  the  inertial  frame  is: 


T  r  V  t  s^blN 

JJtO  +(0  Xjfc(0 


(33) 


As  expected,  leaving  all  terms  in  the  body  frame  and  taking  the  cross  product 

results  in  the  original  equations  for  a  rigid  body  system. 
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4,  =  [*4>y]  +  ~  J b-,  4 

4,  =  [y®2  ]  +  (*4n  —  y )  ^*3  Y  (34) 

43  +  (°42  —  "4ii)®l®2  ^3 

D.  EQUATIONS  OF  MOTION  FOR  THE  GIMBAE 

As  defined  in  Equation  (27),  the  angular  momentum  of  the  gimbal  is: 

hg=Jg«gW  (35) 

In  this  case,  both  J  and  or  '  are  defined  in  the  gimbal  frame.  Using  the 
transport  theorem  the  inertial  derivative  is: 


dhY  (dhY 

— -  =  — -  +o«"vxh(, 

dt  dt  . 


+  e>gWxJghgMr 


Substituting  0)glN  from  Equation  (28)  in  results  in: 


0)h,N+0)s,b 


+  (ojWJV  +  0)s,b  )  X  Jg  (to*"  +  G)s,b 


)  (37) 


Expanding  the  terms,  Equation  (37)  becomes: 


j(jgc-)Yr^(j^) 


+a)Mxjy/"+o,^xjy/Ar 
+o )hiN  x  jg0y + o)g,b  x  jg(og,b 


For  clarity,  each  of  the  six  terms  will  be  evaluated  separately. 
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1. 


Gimbal  Term  1 


Because  .1  is  constant  in  the  gimbal  frame,  it  can  be  removed  from  the  derivative. 
Then,  the  transport  theorem  is  used  to  shift  back  into  the  body  frame,  which  results  in: 

f  j/t  b/N\\8  f  j(„b!N\\8 

=  J 


dt 


j(jyw)T_  \d{^blN) 


j 


dt 


v 


J 

b/N\~\b 


d(  iOb,N) 


dt 


+  WW8XO)“ 


(39) 


In  order  to  take  the  cross  product  of  0)h/g  and  o)blNwe  recognize  that  ob'8  is 
equivalent  to  -0)s/b ,  which  was  defined  in  Equation  (25)  and  that  both  vectors  need  to  be 

in  the  same  frame  of  reference,  therefore  we  need  to  use  Oih"'N  as  defined  in  Equation 
(3 1).  This  results  in  Equation  (40). 
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2.  Gimbal  Term  2 

Again,  J  is  removed  from  the  derivative  and  because  or /ft  is  defined  in  the 
gimbal  frame  the  result  is  simply: 


dt 


0 

0 


g 


(41) 


3.  Gimbal  Term  3 

As  in  term  1,  we  use  cabg,N  and  take  the  cross  product: 
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4.  Gimbal  Term  4 

We  use  0)slh  andoA  '  in  the  gimbal  frame  and  take  the  cross  product: 


s^g/b  ..  x  ,Jo/N  b  *  x  r^b/N 

CO  xjaco  =  og3  x  J  co 


O' 

J  CO 

ft  ft 
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X 

J„  co„ 

82  82 
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1 

1 _ 

JJ0Js, 


g 


(43) 


5.  Gimbal  Term  5 

We  use  (obg/N  and  as/b  in  the  gimbal  frame  and  take  the  cross  product: 
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(44) 


6.  GimbalTerm6 

Q)glb  is  crossed  with  the  scalar  multiple  of  itself  and  is  therefore  zero: 

<s>8lbY.Sylb=8gly.ig8gl 


"O' 

0 

0 

X 

0 

8 \ 

yj. 

(45) 


These  tenns  can  be  combined  and  grouped  for  the  following  result  for  torque 
resulting  from  the  motion  of  the  gimbal  shown  in  Equation  (46). 
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E.  EQUATIONS  OF  MOTION  FOR  THE  ROTOR 

As  defined  in  Equation  (27),  the  angular  momentum  of  the  rotor  is: 

K=Jr*""  (47) 
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Because  we  assume  the  wheel  is  symmetric  about  rl ,  which  is  collinear  with  gx , 
we  can  define  both  Jr  and  of  in  the  gimbal  frame.  Again,  using  the  transport  theorem 
the  inertial  derivative  is: 


fdh ) 

N 

rdh  \ 

— 

r 

v  dt  , 

^  dt  y 

d(  Jf0)''/iV) 


+  (o''/N  x  h(. 

r/N\\g 


dt 


(48) 


I  T  ,-r/JV 

+  0)  xjir« 


J 


Substituting  t orlN  from  Equation  (28)  includes  the  following  nine  terms: 


N 

f d(ir\ 

(toWiV +(Og/A 

\  dt  j 

V 

dt 

)) 


+  UWJV+o>g/i  +  ar/* 


H( 


(Ob,N  +(Og/A  +  G>r/* 


)  (49) 


Expanding  the  terms,  Equation  (37)  becomes: 


dt 


rfK")Y , ,  fdK)T , .  fdK‘) 


+  J, 


dt 


+  J 


dt 


+  ((e>4W  x  JrmblN )  +  (o)g/6  x  Jr©4/JV )) 

+  ((e>4/JV  x  Jr©*/6 )  +  ((og/b  x  Jr©*/4 )) 

+  ((e>4/*  x  Jror/g )  +  (ra*/4  x  Jro)r,g )) 


(50) 


With  the  exception  of  the  change  in  the  inertia  dyad,  Rotor  Terms  1,  2,  4,  5,  6, 
and  7  are  identical  to  Gimbal  Terms  1  through  6.  By  the  same  evaluation  as  shown  in 
section  III.D,  Equations  (51)  through  (56)  follow. 


1. 


Rotor  Term  1 


dt 
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(51) 
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2. 


Rotor  Term  2 


c  J(wg/6)Y 

dt 
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Rotor  Term  4 
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Rotor  Term  5 


0)s'b  x  J  g>0/jv  = 
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Rotor  Term  6 
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6.  Rotor  Term  7 

«g/6xjge)^=0 

Rotor  Term  3,  8  and  9  will  be  evaluated  in  full. 


(52) 


(53) 


(54) 


(55) 


(56) 


7.  Rotor  Term  3 

Again,  J(.  is  removed  from  the  derivative  and  because  e)'/gis  defined  in  the 
gimbal  frame  the  result  is  simply: 
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Rotor  Term  8 


Using  d)  ’"is  in  the  gimbal  frame  and  cbr'g  ,  the  cross  product  becomes: 
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9.  Rotor  Term  9 

Both  (os'b  and  to'/g  are  defined  in  the  gimbal  frame  so  the  cross  product  is: 


oU/AxJ(.o)'Vg  = 
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Using  the  previously  stated  assumption  that  the  rotor  is  symmetric  about  the 
gj  axis,  we  assume  that  J  equals  Jr  ,  and  the  simplified  equation  becomes  Equation  (60) 


26 


Jr® 1 

J>Mi  1 

0  " 

~jr  n" 

r\ 

K  = 

J, A_ 

b  + 

~JrMt 

0 

g  + 

0 

g  + 

0 

0 

g 


0 

-J  8(0 

J  r  8(o  ' 

>l  Si 

r3  Si 

+ 

<v\(7-4) 

g  + 

J,  8(0 

r\  &l 

g  + 

-J,  8(o 

r3  8 1 

0 

0 

g 


0 

.7  Cl  co. 


Si 


-J.Cla ) 


S  1 


1 

O 

i _ 

g  + 

Jfij 

r\ 

i 

o 

_ 1 

g 


(60) 


F.  EQUATIONS  OF  MOTION  FOR  THE  SYSTEM 

Substituting  Equations  (34),  (46)  and  (60)  into  Equation  (24)  the  full  equations  of 
motion  for  a  single  CMG  attached  to  a  rigid  body  can  be  written  as  Equation  (61). 
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In  the  original  equation  fonn,  the  reduced  model  referenced  earlier  in  this 
document  is  shown  in  Equation  (62). 
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t cmg  =  (4  +  4  )  <?£i  +(wW"  x  hCMG  )  +  («g/*  x  hCMG  )  (62) 

where  hCMG  =  Jr0)r/g .  Based  on  Equation  (61),  the  reduced  model  can  also  be 


written  as  Equation  (63). 
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We  make  the  notational  simplification  that: 
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Clearly,  this  equation  could  be  simplified;  however,  because  the  purpose  of 
modeling  is  to  determine  the  relative  significance  of  the  various  terms,  it  will  be  left  in  its 
expanded  form.  This  equation  remains  the  effect  of  one  CMG  only.  In  order  to  increase 
the  fidelity  of  the  model  and  obtain  results  that  are  comparable  to  experimental  data,  it  is 
necessary  to  expand  the  equations  to  an  N-CMG  constellation.  A  4-CMG  system  is 
shown  in  Figure  6  as  an  example  of  a  typical  control  system. 


Figure  6.  Four  CMG  pyramid  configuration.  Both  </>  (offset)  and  /?  (skew)  angles  are 
fixed,  determined  by  the  mounting  of  the  CMGs  with  respect  to  the  spacecraft 

reference  frame  ( bx ,  b2 ,  b3 ). 

While  the  majority  of  the  CMG  terms  have  been  left  in  the  gimbal  frame  of 
reference,  the  transformation  matrices  are  required  both  to  convert  the  angular  velocity  of 
the  body  into  the  gimbal  frame  and  to  convert  the  gimbal  torque  into  the  body  frame  for 
integration. 

Figure  6  shows  that  each  CMG  undergoes  three  rotations  between  the  gimbal 
frame  of  reference  and  the  body  frame  of  reference.  Initially  each  CMG  axis  is  aligned 
with  the  body  axis,  therefore  the  first  rotation  is  a  rotation  about  the  g3  axis  to  create  an 
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offset  angle  between  each  CMG.  In  this  model,  CMG#1  rotates  90°,  CMG#2  rotation 
180°,  CMG#3  rotates  270°  and  CMG#4  rotates  360°.  Then  each  CMG  is  rotated  about 
its  spin  axis  (gj)  an  angle  of  p.  This  is  the  baseline  position  of  each  CMG.  In  order  to 

create  torque,  each  CMG  then  rotates  an  angle  8  about  its  own  g3  axis.  The  rotation 
matrices  are  therefore: 

cos(^)  -sin(^)  0 
DCM,,  (*)  =  sin(^)  cos(^)  0 
0  0  1 

"l  0  0 

dcm4,(/?)=  0  cos  (/?)  -sin(/?)  I  (66) 

0  sin(/?)  cos(/?) 

cos(J)  -sin(^)  0 
DCM^(£)  =  sin(^)  cos(c7)  0 
0  0  1 

The  complete  transfonnation  matrix  from  the  body  frame  to  the  gimbal  frame  is 

DCM,,  W-DCM,,  (A>DCMfc  (*) 

Or 

c8c(j)-%8c(:h(j>  cSs(j>+sSc/3c<j>  sSs/3 
CGIB  =  -sSc<f>-cSc/3s<f>  -%8%(/>+c8c(3c(j)  ec5s/? 

s/3s</>  -s  f3c(f)  c/3 

where  =  cos  (£),.$£  =  sin  (8 ),  c</>  =  cos  (</>),  s</>  =  sin  (</>),  c/3  =  cos(/?), 

s/3  =  sin (/? ) .  Using  the  relationship: 

[c^]"  =  [cs/"]r  =  c^ 

the  DCM  from  the  gimbal  frame  to  the  body  frame  is: 

c8c</>-s8c/3s(/>  -sSctp-cScpstj)  s/3s</> 

CBlG  =  cSs</>+sSc/3c<f>  -sSs<f>+c8c/3c</>  -s/3c</> 
s8s/3  c8s/3  c/3 


(67) 


(68) 

and 


(69) 


(70) 
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27  Bong  Wie,  Space  Vehicle  Dynamics  and  Control,  2nd  ed.  (Reston,  VA:  American  Institute  of 
Aeronautics  and  Astronautics,  2008),  332,  403. 
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IV.  VERIFICATION  OF  MATHEMATICAL  MODELS 


A.  SIMULINK  MODEL 

In  order  to  verify  the  equations,  a  SIMULINK  model  was  created,  which 
evaluated  the  equations  given  control  inputs  o  IV)' .  A  function  block  was  created  that 
accepted  as  inputs  the  body  quaternions  ( q  )  and  angular  velocities  ( (o )  and  the  CMG 

input  of  8  The  block  was  programmed  to  output  both  <7  and  6)  as  well  as  a  variety  of 
intermediate  values  for  CMG  torque  contributions  and  angular  momentums  of  the 
system.  The  function  block  accomplished  the  various  frame  of  reference  transformations 
and  calculated  each  tenn  of  Equations  (61)  and  (65)  for  further  analysis.  Shown  in 
Figure  7,  it  was  used  to  accomplish  both  open  and  closed  loop  simulations. 
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Figure  7.  System  Dynamics  function  block. 
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The  open  loop  simulations  involved  directly  inputting  values  into  the  function 
block  in  order  to  calculate  the  CMG  response. 


The  closed  loop  simulation  used  an  Eigenaxis  Rotational  Maneuver  controller  to 
generate  a  required  control  torque  in  the  body  frame.  The  quaternion  error  was 
calculated  using  Equation  (71)28 
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where  qc  is  the  commanded  quaternion  vector  and  q  is  the  body  quaternion.  The  control 
torque  was  then  generated  using  Equation  (72).29 

u  =  -Kp I;,qe  [1  2  3]  -  KdIb(o  +  a  x  Ib(o  (72) 

A'/  and  Kd  are  the  proportional  and  derivative  gains,  respectively,  calculated  using 
Equations  (73)  and  (74)  and  qe  [l  2  3]  indicates  that  only  the  first  three  elements  of  q 
are  used  in  the  calculations. 

Kp=co2n  (73) 


Kd=2conC 


(74) 


The  control  torque  was  then  used  to  calculate  the  required  gimbal  angle,  rate  and 
acceleration.  Gimbal  rate  was  calculated  using  pseudo-inverse  steering  logic  with  null 
motion  singularity  avoidance. 


as 


From  Equations  (23)  and  (24),  the  external  torque  on  the  spacecraft  can  be  written 


<*(»„))"  _G(h,„1) 


dt 


dt 


+  0)“xH, 


(75) 


The  total  angular  momentum  from  Equation  (23)  can  be  restated  as 


H, 


Hb+hCMG  Jb***  ^^CMG 


(76) 


where  hCMG  =  H  +  Hr .  Combining  Equations  (75)  and  (76)  we  obtain 
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text  ~  (' Jb“’  N  +  ^CMG  )  +  V  X  (' Jb^’  V  +  hCMG  ))  (77) 

The  control  torque  uis  defined  as  the  torque  generated  by  the  CMGs  and  is 
therefore 

u  — — (hcMG+t®  x^cmg)  (7^) 

For  the  purposes  of  developing  an  equation  for  gimbal  rate,  the  angular 
momentum  of  the  CMG  system  is  simplified  beyond  even  the  reduced  model.  The 
complete  angular  momentum  of  the  CMG  is  assumed  to  be 

^cmg  =  ^r^Si  (79) 

For  a  system  of  n-CMGs,  and  defining  h0  =  J  Q  the  angular  momentum  of  the  CMGs  in 
the  body  frame  would  be 

cSpft-sSpfisfy 

hcMG^  =  Y,K  \  I  (80) 

i=1  [  sStsj3 

The  transformation  vector  is  the  first  column  of  the  DCM  from  the  gimbal  frame 
to  the  body  frame  ( CB/G  ).  Having  defined  the  angular  momentum  of  the  CMG  thus,  the 
time  derivative  of  the  angular  momentum  vector  is  therefore 

—sSjc<f>j  -  cScfh^ 

hcMG^Z^o  4  -sd^+cdppcfr  (81) 

,=i  L  cS>sP 

where  the  transformation  vector  is  the  derivative  of  the  first  column  of  the  DCM  from  the 
gimbal  frame  to  the  body  frame.  If  we  define 

-sStc t/i  -  c8)cfh(f)) 

Uj  =  -sSfify  +C  c>.c  (82) 

cS^/3 

and 

A  =  /?0[ai|a2|---an]  (83) 

Equation(81)  becomes 
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h  =  A5 


(84) 


Combining  Equation(78)  and  (84)  the  gimbal  rate  can  be  obtained  as 

8  =  A"1(-u-©6/JvxhCMG)  (85) 

which  was  used  in  Section  I  as  Equation  (11).  Equation  (85)  is  only  valid  for  n  =  3, 
therefore  a  more  general  approach  uses  the  pseudo-inverse  shown  in  Equation  (86). 

8  =  (Al’(AAr)  1)(-u-«6/wxhCMG)  (86) 

A  known  issue  with  obtaining  the  gimbal  rates  using  a  pseudo-inverse  technique 
such  as  Equation  (86)  is  the  existence  of  singular  states.  A  singular  state  “occurs  when 
all  the  individual  CMG  torque  output  vectors  a.  are  perpendicular  to  the  commanded 

torque  direction.”30  One  technique  to  avoid  singularities  is  to  use  null  motion  added  to 
Equation  (86)  such  that  the  gimbal  rates  becomes 

8  =  (Ar  ( AAr  +  M)1  )(-u  - (OhlN  x  hCMG )  (87) 

where  the  value  for  A  is  controlled  by  the  singularity  measure.  The  singularity  measure 
is  calculated  using  Equation  (88). 

m  =  ^det(dd7j  (88) 

The  value  for  A  is  calculated  as  shown  in  Equation  (89). 

if  m  >  0.5 

A  =  0 

(89) 

if  m  <  0.5 

A  =  10 

The  commanded  gimbal  rate  was  then  used  to  calculate  commanded  gimbal  angle 
and  acceleration  for  use  in  the  simulation. 

B.  OPEN  LOOP  SIMULATIONS 

Several  open  loop  Monte  Carlo  simulations  were  conducted  where  random  values 
were  generated  for5,  to,  Q,  Jr,  and  J,r .  The  values  were  bounded  as  shown  in  Table  1. 
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Maximum  Value 

Minimum  Value 

8 

10  rad/ sec" 

2 

-10  rad/sec 

CO 

0.5  rad/sec  (-28.5  °/sec) 

-0.5  rad/sec  (-28.5  °/sec) 

Q 

1000  rpm 

15,000  rpm 

Jr 

HW 

rotor  scaling 

"1.25  0  0 

0  0.625  0 

0  0  0.625 

kg-m2 

HW 

rotor  scaling 

“0.75  0  0 

0  0.375  0 

0  0  0.375 

kg  ■  m2 

Js 

[1.25  0  0  ' 

HW gimbal  scaling  0  1.25  0 

[_  0  0  1.25 

kg  ■m2 

[0.75  0  0 

HWgimhal  scaling  0  0.75  0 

[  0  0  0.75 

kg  ■  m2 

- 

Table  1.  Bounds  for  initial  open  loop  Monte  Carlo  simulation. 


Baseline  moment  of  inertia  tensors  were  chosen  to  simulate  a  representative  CMG 
table  and  a  CubeSat,.  The  system  level  moment  of  inertia  tensor  remained  constant  while 
the  remainder  of  the  inputs  were  bounded  as  shown  in  Table  1.  Inertia  tensor  values  for 
the  CMG  gimbal  and  rotor  were  obtained  from  Honeywell.  Due  to  the  proprietary  nature 
of  the  information,  the  actual  values  are  not  used;  however  normalized  representative 
values  are  used.  The  values  of  the  CMG  inertia  tensors  were  varied  by  25%.  Because 
8  is  the  control  variable  the  bounds  were  set  well  outside  the  typical  values  (10  rad/sec 
equates  to  576  deg/sec  ).  Q  values  were  estimated  using  an  industry  value  for  the  lower 
bound  and  proposed  industry  values  for  the  upper.  The  simulation  was  run  for  one 
million  steps  and  the  output  values  for  the  CMG  gimbal  and  output  torques  were 
calculated  at  each  step. 

A  measure  of  how  accurate  the  reduced  model  was  compared  to  the  full  model 
was  computed  using  the  values  for  the  CMG  Axis  2  (output)  and  Axis  3  (gimbal)  torque 
at  each  step.  The  accuracy  for  each  data  point  was  calculated  using  the  standard  relative 
error,  shown  in  Equation  (90). 
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measure  of  accuracy  = 


(90) 


The  reduced  model  has  no  terms  that  contribute  to  Axis  1 .  Therefore,  the  measure 
of  accuracy  calculation  could  not  be  completed  for  that  accuracy.  However,  in  order  to 
confirm  that  it  is  correct  for  the  reduced  model  to  disregard  Axis  1,  a  measure  of 
significance  was  conducted  where  the  magnitude  of  the  Axis  1  torque  was  compared  to 
the  magnitude  of  the  Axis  2  (output)  torque,  shown  in  Equation  (91). 


measure  of  signifigance  = 


Tfuu  (Axis  1) 

\l_ 

2^  (Axis  2) 

\i 

(91) 


The  accuracy  of  the  entire  model  was  then  estimated  by  taking  the  infinity-norm 
of  the  vectors  generated  by  Equations  (90)  and  (91).  The  infinity  norm  of  the  vector  is 
the  maximum  value  so  the  measure  of  accuracy  in  Table  2  represents  the  worst  case 
scenario  of  the  reduced  model. 


System 

Axis  1 

Measure  of  Significance 

Axis  2 

Measure  of  Accuracy 

Axis  3 

Measure  of  Accuracy 

CubeSat 

0.0138 

0.00208 

0.01195 

CMG  Table 

0.0141 

0.00854 

0.01365 

Table  2.  Open  loop  measure  of  accuracy. 


Table  2  shows  that  Axis  1  torque  is  less  than  1%  of  the  Axis  2  (output)  torque, 
meaning  that  the  reduced  model  failure  to  include  any  magnitude  for  that  Axis  is 
appropriate  for  most  models.  It  also  shows  that  the  reduced  model  is  accurate  to  within 
1%  for  all  input  values.  A  more  detailed  representation  of  the  accuracy  of  the  model  is 
shown  in  Figure  8,  which  shows  the  histograms  of  the  error  vectors  for  the  CMG  table. 
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Number  of  Items  Number  of  Items  Number  of  Items 
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Figure  8.  Histogram  of  open  loop  simulations  for  the  CMG  table.  Note  that  while  there 
are  outliers,  the  model  is  generally  accurate  to  0.00137  for  the  input  torque  and 
0.00086  for  the  output  or  an  order  of  magnitude  better  than  the  worst  outliers. 
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Figure  9  shows  the  histogram  data  for  the  CubeSat  simulations  revealing  similar 


results. 
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Figure  9.  Histogram  of  open  loop  simulations  for  the  CubeSat  values.  Again,  while 
outliers  are  present,  the  model  is  generally  accurate  an  order  of  magnitude  better 

than  the  worst  outliers. 


C.  CLOSED  LOOP  SIMULATIONS 

A  closed  loop  Monte  Carlo  simulation  was  run  where  the  initial  Euler  Angles 
were  varied  from  0-180  and  the  Rotor  and  Gimbal  Inertia  Tensors  were  allowed  to  vary 
by  25%,  similar  to  the  Open  Loop  Monte  Carlo  but  with  Q  values  were  held  constant  and 
the  angular  velocity  and  gimbal  acceleration  detennined  by  the  simulation.  In  each 
simulation  the  commanded  end  state  was  zero  pitch,  roll  and  yaw  and  zero  angular 
velocity.  During  each  simulation,  the  value  of  each  of  the  six  Gimbal  and  nine  Rotor 
terms  were  calculated  for  each  axis  of  each  CMG.  As  shown  in  Figure  10  each  term’s 
magnitude  was  calculated  in  each  CMG  axis. 
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Figure  10.  Value  of  each  of  the  terms  contributing  to  the  overall  CMG  torque. 


At  every  time  step,  the  magnitude  of  each  of  the  terms  was  divided  by  the  tenn 
with  the  maximum  value  to  detennine  the  relative  magnitude  of  each  term  as  shown  in 
Figure  1 1 . 
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Figure  11.  Relative  magnitude  calculation  conducted  at  each  time  step. 


For  each  complete  simulation,  the  maximum  relative  magnitudes  were  saved  for 
comparison.  The  Monte  Carlo  simulation  was  run  for  10,000  iterations.  The  saved  data 
was  analyzed  to  determine  if  there  were  any  conditions  in  which  the  tenns  that  are  not 
included  in  the  reduced  model,  i.e.,  those  terms  believed  to  be  small  enough  to  be 
ignored,  are  significant.  The  maximum  values  for  the  Monte  Carlo  for  the  Gimbal  terms 
and  Rotor  terms  for  Axis  1  torque  are  in  Table  3,  Axis  2  (output)  are  in  Table  4,  and  Axis 
3  (gimbal)  torque  are  in  Table  5. 
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Gimbal  Tenn  1 

Gimbal  Term  2 

Gimbal  Tenn  3 

Gimbal  Term  4 

Gimbal  Tenn  5 

Re( 

g 

el 

bJN  T  bJN 

OJ  *  x  Jo)  * 

oblN  xJgo8/b 

g  dt 

V  7 

dt 

V  7 

luced  Mod 

100.00% 

0.0000% 

100.00% 

100.00% 

100.00% 

Gimbal  Term  6 

Rotor  Term  1 

Rotor  Tenn  2 

Rotor  Term  3 

Rotor  Term  4 

J, 

g 

f^KYY 

f^KYY 

bJN  T  b  IN 

0) 8  x  J  .w  8 

dt 

V  J 

dt 

V  7 

Reduced  Model 

dt 

V  7 

0.0000% 

100.00% 

0.0000% 

0.0000% 

100.00% 

Rotor  Term  5 

Rotor  Tenn  6 

Rotor  Tenn  7 

Rotor  Term  8 

Rotor  Tenn  9 

fi)wxJxras/i 

#/jxjy 

(D^xjy* 

Reduced  Model 

to^xjy* 

Reduced  Model 

58.418% 

62.481% 

0.0000% 

0.0000% 

0.0000% 

Table  3.  Maximum  relative  values  for  Gimbal  and  Rotor  terms  in  Axis  1 . 
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Gimbal  Tenn  1 

Gimbal  Term  2 

Gimbal  Tenn  3 

Gimbal  Term  4 

Gimbal  Tenn  5 

Re< 

g 

el 

bJN  ,  bJN 

(0  s  x  s 

«^xJgw6/" 

0)b,N  X  Jg(Og,b 

dt 

V  7 

dt 

\  7 

iuced  Mod 

1.0527% 

0.0000% 

0.27903% 

1.0197% 

1.0374% 

Gimbal  Term  6 

Rotor  Term  1 

Rotor  Tenn  2 

Rotor  Term  3 

Rotor  Term  4 

G>*/6xJg(o*/6 

J, 

g 

( </(«.*'*  )Y 

pK/g)Y 

bJN  T  bJN 

(8  s  Xjr(0‘ 

dt 

v  J 

dt 

V  7 

0.0000% 

0.29326% 

0.0000% 

0.0000% 

0.27903% 

Rotor  Term  5 

Rotor  Tenn  6 

Rotor  Tenn  7 

Rotor  Term  8 

Rotor  Tenn  9 

(os/bxJr(ob,N 

(0^xjg©^ 

(ob/NxJr(or,g 

Reduced  Model 

(0s,b  x  Jf.tor/g 

Reduced  Model 

0.57009% 

0.30515% 

0.0000% 

100.00% 

100.00% 

Table  4.  Maximum  relative  values  for  Gimbal  and  Rotor  terms  in  Axis  2. 
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Gimbal  Tenn  1 

Gimbal  Term  2 

Gimbal  Tenn  3 

Gimbal  Term  4 

Gimbal  Tenn  5 

Re< 

'd((Og/bp 

g 

el 

bJN  ,  bJN 

CO  x  J^co  5 

(0s/ixJg(0WJV 

wWArxJgw?/6 

dt 

V  J 

dt 

V  ) 

iuced  Mod 

0.0000% 

100.00% 

0.32274% 

0.0000% 

0.0000% 

Gimbal  Term  6 

Rotor  Term  1 

Rotor  Tenn  2 

Rotor  Term  3 

Rotor  Term  4 

G>*/6xJg(0*/6 

J, 

g 

pKYY 

pKYY 

bJN  T  bJN 

CO  xj(0)? 

dt 

\  J 

dt 

V  J 

Reduced  Model 

dt 

V  7 

0.0000% 

0.0000% 

63.145% 

0.0000% 

0.30415% 

Rotor  Term  5 

Rotor  Tenn  6 

Rotor  Tenn  7 

Rotor  Term  8 

Rotor  Tenn  9 

ag/bxJr(oi/N 

(0WNxJg(Os/b 

©g/ixjg(0g/ft 

ob/N  x  Jrco''/g 

Reduced  Model 

(og/b  x  3rorlg 

Reduced  Model 

0.0000% 

0.0000% 

0.0000% 

100.00% 

0.0000% 

Table  5.  Maximum  relative  values  for  Gimbal  and  Rotor  terms  in  Axis  3. 


From  this  information  we  can  see  that  Axis  1  shows  the  most  amount  of  variation 
on  which  terms  contribute  most.  However,  the  open  loop  analysis  showed  that  the 
relative  magnitude  of  the  Axis  1  torque  compared  to  the  Axis  2  torque  makes  the  overall 
output  negligible.  Recall  that  the  reduced  model  assumes  that  the  Axis  1  torque  is  zero. 
For  Axis  3,  the  tenns  that  are  not  included  in  the  reduced  model  do  not  exceed  0.3%  at 
any  point  in  the  simulation.  However,  for  Axis  2,  several  terms  reach  1%  significance 
and  the  accumulation  of  those  tenns  could  potentially  reach  3%.  This  likely  remains 
below  the  threshold  of  significance  for  most  applications  but  could  potentially  be 
meaningful. 
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Section  Notes 


28  Bong  Wie,  Space  Vehicle  Dynamics  and  Control,  2nd  ed.  (Reston,  VA:  American  Institute  of 
Aeronautics  and  Astronautics,  2008),  403. 

29  Ibid.,  406. 

30  Ibid.,  440. 
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V.  EXPERIMENTAL  VERIFICATION 


A.  HONEYWELL  MOMENTUM  CONTROL  SYSTEM 


In  August  201 1,  a  series  of  experimental  maneuvers  was  conducted  on  Honeywell 
Space  System’s  spacecraft  attitude  dynamics  testbed.  The  testbed,  shown  in  Figure  12,  is 
capable  of  high-agility  slewing,  incorporating  six  225  ft-lb  CMGs.31 


Figure  12.  Honeywell’s  Momentum  Control  System  (MCS)  testbed. 

The  testbed  is  capable  of  unlimited  movement  about  the  z-axis,  and  ±30°  in  the  x 

and  y-axes.32  The  vehicle  has  no  umbilical  to  ground,  using  a  wireless  link  for  command 

and  telemetry.33  The  air  bearing  has  a  3200  pound  load  capacity  and  is  mounted  on  a 

concrete  base  that  is  isolated  from  the  surrounding  building.34  For  the  NPS  experimental 

runs,  the  CMG  rotor  speed  was  set  to  1000  rpm.  The  natural  frequency  of  the  testbed, 

used  to  calculate  gains,  was  0.4  Hz  while  the  damping  ratio  was  set  to  0.7. 
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B. 


TRAPEZOID  EXPERIMENTAL  SETUP 


Based  on  input  from  industry  a  standard  trapezoid  maneuver  was  created.  The 
maneuver  consists  of  a  feedforward  body  angular  velocity  and  quaternion  paths.  The 
process  for  calculated  an  industry  maneuvers  begins  with  establishing  a  desired  set  of 
Euler  Angles.  The  example  maneuver  used  was  from  +2.25°  pitch,  —20.25°  roll  and  0° 
yaw  to  —2.25°  pitch,  +20.25°  roll  and  0°  yaw.  The  desired  Euler  Angles  along  with  the 
maximum  body  angular  velocity  and  maximum  body  acceleration  are  used  to  generate 
angular  velocity  paths  that  resemble  trapezoids,  as  shown  in  Figure  13. 


Commanded  Angular  Velocity  (oo) 


Figure  13.  Example  of  an  industry  standard  trapezoid  style  maneuver. 

The  Commanded  Euler  Angles  were  also  used  to  generate  Commanded 
Quaternion  paths,  as  shown  in  Figure  14. 
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Commanded  Quaternions 


Figure  14.  Commanded  quaternions  for  an  industry  standard  maneuver. 

C.  COMPARISON  OF  RESULTS 

The  results  of  the  simulation  were  compared  in  several  ways.  The  simulation 
results  were  compared  to  the  commanded  angular  velocity  and  quaternions  to  determine 
how  well  the  simulation  tracked  the  feedforward  commands  as  shown  in  Figure  15  and 
Figure  16. 
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Figure  15.  Comparison  of  commanded  angular  velocity  to  simulation  output  angular 

velocity. 
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State  Trajectories  (q) 


time  (seconds) 


Figure  16.  Comparison  of  commanded  quaternions  to  simulation  output  quaternions. 

These  figures  show  that  the  simulation  tracks  the  input  extremely  well.  This  was 
confirmed  by  calculating  the  Root  Mean  Square(RMS)  error  of  the  commanded  input 
versus  the  output.  The  RMS  error  for  the  angular  velocity  is  1. 9268  xlCT5  and 
5.3755 xlO  6  for  the  quaternions. 

The  simulation  results  were  also  compared  to  the  Honeywell  MCS  output  of  the 
same  commanded  angular  velocity  and  quaternions  to  determine  how  well  the  simulation 
reflects  the  hardware  as  shown  in  Figure  17  and  Figure  18. 
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Figure  17.  Comparison  of  simulation  output  angular  velocity  to 

hardware  output  angular  velocity. 


State  Trajectories  (q) 
(hardware  to  simulation  comparison) 


time  (seconds) 

Figure  18.  Comparison  of  simulation  output  quaternions  to 

hardware  output  quaternions. 


Figure  17  shows  that  the  physical  hardware  has  transient  responses  not  revealed  by  the 
simulation.  These  reflect  inaccuracies  in  the  transient  response  of  the  simulation.  For 
instance;  the  simulation  does  not  currently  account  for  friction  effects  on  the  individual 
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CMG  gimbals.  Even  with  this  transient  error,  the  overall  RMS  error  remains  low.  The 
RMS  error  for  the  angular  velocity  is  5.6897 xlCT4;  for  the  quaternions  it  is 
9.0234X1CT4. 

The  Honeywell  MCS  data  includes  most  of  the  information  needed  to  calculate 
the  values  of  the  individual  terms.  The  gimbal  acceleration  is  not  captured  directly  but 
gimbal  torque  is  measured.  Since  all  the  other  terms  can  be  calculated,  we  can  use  the 
gimbal  torque  to  solve  for  Gimbal  Term  2  and  Rotor  Term  2.  The  maximum  relative 
magnitudes  of  each  tenn  were  calculated  over  the  course  of  the  simulation  in  the  same 
manner  as  described  in  the  closed  loop  analysis.  Table  6  contains  the  maximum  relative 
values  for  Axis  1 .  While  this  axis  shows  the  most  variation,  the  total  magnitude  of  the 
torque  in  this  axis  remains  small  in  relation  to  the  Axis  2  or  output  torque. 


Gimbal  Tenn  1 

Gimbal  Term  2 

Gimbal  Tenn  3 

Gimbal  Term  4 

Gimbal  Term  5 

x 

Re< 

g 

el 

bJN  ,  bJN 

(0  8  x  J^ra  8 

<n*lbxJg&b/N 

w“xjg0)?/6 

dt 

V  7 

90.5840% 

0.0000% 

100.0000% 

72.6682% 

100.0000% 

Gimbal  Term  6 

Rotor  Term  1 

Rotor  Tenn  2 

Rotor  Term  3 

Rotor  Term  4 

G>*/6xJg(0*/6 

J, 

g 

,  GKOY 

t  (d(^)  y 

bJN  T  h/N 

Cl)  xJo8 

dt 

V  J 

dt 

V  7 

0.0000% 

54.2066% 

0.0000% 

0.0000% 

100.0000% 

Rotor  Term  5 

Rotor  Tenn  6 

Rotor  Tenn  7 

Rotor  Term  8 

Rotor  Tenn  9 

rag/fe  x  J;V/W 

W6WxJgO)g/6 

(Og/fe  x  Jg<Og,b 

0)iWxJ,,0)r/? 

Reduced  Model 

e>s/bxj  rtor,g 

Reduced  Model 

36.0455% 

27.2439% 

0.0000% 

0.0000% 

0.0000% 

Table  6.  Maximum  relative  values  for  Gimbal  and  Rotor  terms  in  Axis  1 . 


Table  7  shows  the  results  for  Axis  2.  As  expected  the  terms  not  included  in  the 
reduced  model  are  negligible. 
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Gimbal  Tenn  1 

Gimbal  Term  2 

Gimbal  Tenn  3 

Gimbal  Term  4 

Gimbal  Term  5 

Je 

Re< 

r  d  ^ 

g 

el 

bJN  ,  bJN 

CO  x  J^co  s 

(0s/ixJg(0WJV 

(Ob/NxJg(Og/b 

dt 

V  J 

dt 

V  ) 

iuced  Mod 

0.0059% 

0.0000% 

0.0000% 

0.0059% 

0.0059% 

Gimbal  Term  6 

Rotor  Term  1 

Rotor  Tenn  2 

Rotor  Term  3 

Rotor  Term  4 

G>*/4xJg(0*/4 

J, 

g 

(d(i »^)Y 

bJN  T  bJN 

(0  xj(0)? 

dt 

\  J 

dt 

V  J 

Reduced  Model 

dt 

V  J 

0.0000% 

0.0017% 

0.0000% 

0.0000% 

0.0017% 

Rotor  Term  5 

Rotor  Tenn  6 

Rotor  Tenn  7 

Rotor  Term  8 

Rotor  Tenn  9 

ro^xjy^ 

0)b/N  x  Jg(Os/b 

©g/ixjg(0g/ft 

ob/N  xJrar/8 

Reduced  Model 

(o8/b  x  Jf.tor/g 

Reduced  Model 

0.0034% 

0.0017% 

0.0000% 

100.0000% 

100.0000% 

Table  7.  Maximum  relative  values  for  Gimbal  and  Rotor  terms  in  Axis  2. 


Table  8  shows  the  results  for  Axis  3.  Again,  the  terms  not  included  in  the  reduced 
model  are  negligible. 
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Gimbal  Tenn  1 

Gimbal  Term  2 

Gimbal  Tenn  3 

Gimbal  Term  4 

Gimbal  Term  5 

Je 

Re< 

r  d  ^ 

g 

el 

bJN  ,  bJN 

CO  x  J^co  s 

(0s/ixJg(0WJV 

(0b/NxJg(08/b 

dt 

V  J 

dt 

V  ) 

iuced  Mod 

0.0000% 

100.0000% 

0.0000% 

0.0000% 

0.0000% 

Gimbal  Term  6 

Rotor  Term  1 

Rotor  Tenn  2 

Rotor  Term  3 

Rotor  Term  4 

G>*/4xJg(0*/4 

J, 

g 

(d(i »^)Y 

bJN  T  bJN 

(0  xj(0)? 

dt 

v  J 

dt 

V  J 

Reduced  Model 

dt 

V  J 

0.0000% 

0.0000% 

38.0356% 

0.0000% 

0.0017% 

Rotor  Term  5 

Rotor  Tenn  6 

Rotor  Tenn  7 

Rotor  Term  8 

Rotor  Tenn  9 

ro^xjy^ 

0)b/N  x  Jg(Os/b 

©g/ixjg(0g/ft 

ob/N  xJrar/8 

Reduced  Model 

(o8/b  x  Jf.tor/g 

Reduced  Model 

0.0000% 

0.0000% 

0.0000% 

100.0000% 

0.0000% 

Table  8.  Maximum  relative  values  for  Gimbal  and  Rotor  terms  in  Axis  3. 

D.  OPTIMAL  MANEUVER  EXPERIMENTAL  SETUP 


In  addition  to  the  industry  standard  maneuver,  the  same  initial  and  end  Euler 
Angles  were  used  to  develop  a  time  optimal  maneuver.  The  time  optimal  maneuver 
commanded  angular  velocity  is  shown  in  Figure  19  and  the  commanded  Quaternions  are 
shown  in  Figure  20. 
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Angular  Velocity  (rad/sec)  Angular  Velocity  (rad/sec) 


Commanded  Angular  Velocity  (co) 


Figure  19.  Time  optimal  commanded  angular  velocity. 


Commanded  Quaternions 


E.  COMPARISON  OF  RESULTS 

The  results  of  the  simulation  were  compared  in  the  same  manner  as  the  trapezoid 
maneuver.  The  simulation  results  were  compared  to  the  commanded  angular  velocity 
and  quaternions  to  determine  how  well  the  simulation  tracked  the  feedforward  commands 
as  shown  in  Figure  21  and  Figure  22. 


55 


State  Trajectories  (to) 


Figure  2 1 .  Comparison  of  commanded  angular  velocity  to 
simulation  output  angular  velocity. 


State  Trajectories  (q) 

(input  to  simulation  comparision) 


time  (seconds) 


Figure  22.  Comparison  of  commanded  quaternions  to 
simulation  output  quaternions. 


These  figures  show  differences  between  the  simulation  input  and  output  are 
negligible.  This  was  again  confirmed  by  calculating  the  RMS  error  of  the  commanded 
input  versus  the  output.  The  RMS  error  for  the  angular  velocity  is  1. 3003  xlO  3  and 
6.55 17  x  1CT4  for  the  quaternions. 

The  simulation  results  were  also  compared  to  the  Honeywell  MCS  output  of  the 
same  commanded  angular  velocity  and  quaternions  to  determine  how  well  the  simulation 
reflects  the  hardware  as  shown  in  Figure  23  and  Figure  24. 
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State  Trajectories  (oo) 


Figure  23.  Comparison  of  simulation  output  angular  velocity  to  hardware  output  angular 

velocity. 


State  Trajectories  (q) 
(hardware  to  simulation  comparison) 


Figure  24.  Comparison  of  simulation  output  quaternions  to  hardware  output  quaternions. 


Figure  23Figure  23.  again  shows  transient  responses  in  the  hardware  not  revealed  by  the 
simulation.  The  transient  differences  are  greater  in  the  optimal  control  maneuver  than  in 
the  trapezoid  maneuver  due  to  the  increased  aggressiveness  of  the  optimal  control  path. 
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This  is  also  reflected  in  increased  RMS  error.  The  RMS  error  for  the  angular  velocity  is 
3.5693xlCT3  and  1 .9205 xlCT3  for  the  quaternions. 

The  maximum  relative  magnitudes  of  each  term  were  calculated  over  the  course 
of  the  simulation  in  the  same  manner  as  the  trapezoid  maneuver.  Table  9,  Table  10,  and 
Table  11  contain  the  maximum  relative  values  for  Axis  1,  2  and  3,  respectively.  The 
computed  values  are  consistent  with  the  expected  results  and  with  the  trapezoid 
maneuver. 
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Table  9.  Maximum  relative  values  for  Gimbal  and  Rotor  terms  in  Axis  1. 


58 


Gimbal  Tenn  1 

Gimbal  Term  2 

Gimbal  Tenn  3 

Gimbal  Term  4 

Gimbal  Term  5 

t  fd(^N))8 

Re< 

"</(tog/i)^ 

g 

el 

bJN  ,  bJN 

(0  X  J^tD  8 

<os/b  X  Jgto6W 

ab/NxJgag/b 

Jg  dt 

V  J 

dt 

v  9 

iuced  Mod 

0.0280% 

0.0000% 

0.0000% 

0.0280% 

0.0280% 

Gimbal  Term  6 

Rotor  Term  1 

Rotor  Tenn  2 

Rotor  Term  3 

Rotor  Term  4 

G>*/6xJg©*/6 

J, 

'd(e>b/"p 

g 

( </(«."*  )Y 

p(^/g)T 

bJN  T  bJN 

to 8  xJ  co8 

dt 

\  J 

dt 

V  J 

0.0000% 

0.0081% 

0.0000% 

0.0000% 

0.0068% 

Rotor  Term  5 

Rotor  Tenn  6 

Rotor  Tenn  7 

Rotor  Term  8 

Rotor  Tenn  9 

ro^xjy^ 

0)b,N  X  Jg0)s,b 

tog/*  X  Jg(Og,b 

oblN  x  Jr(or/g 

Reduced  Model 

as,b  x  Jrt»r/* 

Reduced  Model 

0.0162% 

0.0081% 

0.0000% 

100.0000% 

100.0000% 

Table  10.  Maximum  relative  values  for  Gimbal  and  Rotor  terms  in  Axis  2. 
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Table  1 1 .  Maximum  relative  values  for  Gimbal  and  Rotor  terms  in  Axis  3. 


F.  NPS  SATELLITE  SIMULATOR  TEST  BED  (STB) 

The  NPS  Satellite  Simulator  Test  Bed  (STB)  shown  in  Figure  25  is  an  Andrews 
Space  Satellite  Simulator  tailored  for  use  at  NPS.  Four  seven  ft-lb  CMGs  are  mounted  in 
a  roof-down  configuration  on  an  attitude  platform,  which  is  supported  by  a  spherical  air 
bearing. 
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Figure  25.  NPS  Satellite  Simulator  Test  Bed  (STB). 


The  simulator  is  capable  of  unlimited  movement  about  the  z-axis,  and  +45°  in  the 
x  and  y-axes.  The  rotor  speed  was  set  to  5000  rpm  for  the  experimental  maneuver. 

G.  DESCRIPTIONS  OF  COMPARISON  MANEUVER 

Currently  the  STB  is  commanded  using  feedforward  quaternions  only. 
Additionally,  CMG1  has  malfunctioned  and  has  been  shut  down,  leaving  only  three 
functional  CMGs.  The  test  maneuver  commanded  quaternions  are  shown  in  Figure  26. 
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Commanded  Quaternions 


Figure  26.  STB  commanded  quaternions. 

H.  COMPARISON  OF  RESULTS 

Again,  the  results  of  the  simulation  were  compared  in  several  ways.  The 
simulation  results  were  compared  to  the  commanded  quaternions  to  determine  how  well 
the  simulation  tracked  the  feedforward  commands  as  shown  in  Figure  27. 


State  Trajectories  (q) 

(input  to  simulation  comparision) 


Figure  27.  Comparison  of  commanded  quaternions  to  simulation  output  quaternions. 

Figure  27  shows  that  the  simulation  exhibits  a  slight  delay  before  tracking  the 
commanded  quaternions.  The  delay  may  be  caused  by  a  lack  of  expected  feedforward 
angular  velocity.  This  also  leads  to  an  increased  RMS  error  for  this  data  of  1. 1970  x  1CT2 . 


62 


The  simulation  results  were  also  compared  to  the  STB  output  of  the  same 
commanded  quaternions  to  determine  how  well  the  simulation  reflects  the  hardware  as 
shown  in  Figure  28. 


State  Trajectories  (q) 


Figure  28.  Comparison  of  simulation  output  quaternions  to  hardware  output  quaternions. 

Figure  28  reveals  that  the  hardware  response  lags  slightly  behind  the  simulation. 
The  causes  of  this  lag  were  not  investigated  but  likely  reflect  the  imperfections  in  the 
model  discussed  above.  This  is  reflected  in  an  increased  RMS  error  of  2.8 193  x  10  2 . 

The  STB  data  includes  most  of  the  information  needed  to  calculate  the  values  of 
the  individual  terms.  However,  it  is  not  instrumented  to  capture  either  gimbal 
acceleration  or  gimbal  torque.  This  lack  of  data  prevents  calculation  of  Gimbal  Term  2 
and  Rotor  Term  3,  however,  these  terms  are  assumed  to  be  significant,  therefore  the 
inability  to  determine  a  relative  magnitude  is  acceptable.  The  maximum  relative 
magnitudes  of  each  term  were  again  calculated  over  the  course  of  the  simulation  in  the 
same  manner  as  the  closed  loop  analysis  above. 
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Table  12.  Maximum  relative  values  for  Gimbal  and  Rotor  terms  in  Axis  1. 
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Table  13.  Maximum  relative  values  for  Gimbal  and  Rotor  terms  in  Axis  2. 
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Table  14.  Maximum  relative  values  for  Gimbal  and  Rotor  terms  in  Axis  3. 


Table  12,  Table  13,  and  Table  14  confirm  the  results  of  the  previous  experimental 
analysis. 
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Section  Notes 


3 1  Mason  A.  Peck  and  Susan  C.  Kim,  "New  Results  from  the  Spacecraft  Momentum  Control  and 
Line-of-Sight  Testbed,"  Advances  in  the  Astronautical  Sciences  118  (2004),  275. 

32  Brian  Underhill  and  Brian  Hamilton,  "Momentum  Control  System  and  Line-of-Sight  Testbed 
(AAS  06-053)"  (Advances  in  the  Astronautical  Sciences  Guidance  and  Control,  Breckenridge,  CO, 
February  4-8,  2006),  554. 

33  Ibid.,  554. 

34  Ibid.,  554 
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VI.  CONCLUSIONS  AND  FUTURE  WORK 


A.  CONCLUSIONS 

The  high  fidelity  model  of  a  CMG  controlled  spacecraft  results  in  a  series  of 
equations  that  significantly  increase  the  complexity  compared  to  current  control  models. 
The  results  from  the  verification  and  the  experimental  analysis  show  that  the  reduced 
model  is  highly  accurate  and  can  be  successfully  used  to  generate  optimal  path  solutions. 
Specifically,  the  reduced  model  is  accurate  to  greater  than  0.01%  in  the  gimbal  torque 
axis.  The  closed  loop  analysis  also  reveals  that  several  terms  in  the  output  axis  may 
reach  significance  in  certain  cases.  During  closed  loop  testing  the  coupling  tenns 
between  the  angular  velocity  of  the  CMG  and  the  angular  velocity  of  the  body 
(tog/i  x  Jgwb/N  and  to*  v  x  J_of  reach  1%  relative  value.  With  the  addition  of  other, 

less  significant  terms,  the  total  could  reach  3%  of  the  overall  output  torque.  In  the  effort 
to  maximize  the  ability  to  slew  quickly  and  accurately,  these  terms  may  be  exploitable. 

B.  FUTURE  WORK 

The  current  model  is  accurate  in  the  steady  state  but  inaccurate  in  the  transient, 
best  illustrated  in  Figure  17,  and  could  be  improved  by  including  more  aspects  of  the 
physical  system,  including  friction  terms  and  testbed  imbalances.  Specifically,  the  NPS 
testbed  has  an  offset  between  its  center  of  gravity  and  center  of  rotation.  Including  this 
offset  in  the  model  would  increase  the  accuracy  of  the  output. 

While  a  large  amount  telemetry  data  is  available  from  the  Honeywell  MCS,  a 
more  instrumented  CMG  /  CMG  table  would  increase  the  ability  to  determine  when 
CMG  back  drive  becomes  significant  and  how  it  can  avoided.  This  would  also  assist 
with  increasing  the  accuracy  of  the  model.  All  of  this  work  could  then  be  used  to 
investigate  how  the  high  fidelity  model  could  be  used  to  calculate  increased  optimal 
solutions  to  the  satellite  slewing  problem. 
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